#########################################################################################
# Title: Figure5.R
# Author: Yuki Atsusaka (atsusaka@rice.edu)
# Aim: Code to replicate Figure 5 in Atsusaka (2021)
#########################################################################################

# CLEAN THE GLOBAL ENVIRONMENT AND READ PACKAGES
rm(list=ls())
library(scales)

################################################################
# DEFINING FUNCTIONS TO DRAW PANELS A-D IN FIGURE 5
################################################################

#++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++#
# FUNCTION TO PREDICT THE PROBABILITY OF MINORITY CANDIDATE EMERGENCE AND VICTORY
#++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++#

minorep <- function(M, C){

q = sqrt(M*C) - 50 # GEOMETRIC MEAN
p = pnorm(q=q,  mean=0, sd=1)
  
return(p)
}

#++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++#
# FUNCTION TO SIMULATE THE RACIAL MARGIN OF VICTORY
#++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++#

M.sim <- function(C, coethnic, crossover){
   
V_m = C*coethnic + ((100-C)*crossover)
V_w = C*(1-coethnic) + ((100-C)*(1-crossover))
M = 1/2*(V_m - V_w) + 50

return(M)   
}


#++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++#
# FUNCTION TO PREDICT THE IMPACT OF REDISTRICTING
#++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++#

redistrict <- function(coethnic, crossover, gap=NULL){

# PERCENTAGE OF MINORITY VOTERS
C = seq(from=1, to=100, by=0.1) 
    
# TURNOUT ADJUSTED PERCENTAGE OF MINORITY VOTERS
if(is.null(gap)){
C = C
}else{
C = (C*gap[1])/(C*gap[1] + (100-C)*gap[2])*100
}

# SIMULATED RACIAL MARGIN OF VICTORY
V_m = C*coethnic + ((100-C)*crossover)
V_w = C*(1-coethnic) + ((100-C)*(1-crossover))
M = 1/2*(V_m - V_w) + 50 

# MODEL PREDICTIONS
MC.sqrt.min50 = sqrt(C*M) - 50  # Geometric Mean 
out <- pnorm(MC.sqrt.min50)     # Probability of Candidate Emergence

return(out)
}


#++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++#
# FUNCTION TO SPREDICT THE NUMBER OF MINORITY WINNERS
#++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++#

n.minorep <- function(out){

# MONTE CARLO SIMULATION
N_pred <- NA

for(i in 1:1000){
model.sample = sapply(out, function(x) rbinom(n=1, size=1, prob=x))
N_pred[i] <- sum(model.sample)  
} # CLOSE THE LOOP

return(N_pred)
}

   
############################################################ 
# CREATING AND COMBINING FOUR PANELS IN FIGURE 5
############################################################ 


pdf(here::here("Atsusaka_Figure5.pdf"), height=8, width=8) # Open a pdf file
par(mfrow=c(2,2), mar=c(4,4,3,1))

#++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++#
# (Panel A) Predicting the Probability with Fixed C over M
#++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++#

C = 40
Mvec <- seq(from=0, to=100, by=0.01) # Adjusted Racial Margin of Victory
Cvec <- rep(x=C, times=length(Mvec))   # Percentage of Minority Voters (SET TO 25)
pred <- minorep(M=Mvec, C=Cvec)

simM1 <- M.sim(C=C, coethnic=1, crossover=0.3)
simM2 <- M.sim(C=C, coethnic=1, crossover=0.4)
simM3 <- M.sim(C=C, coethnic=1, crossover=0.5)

point1 <- minorep(M=simM1, C=C)
point2 <- minorep(M=simM2, C=C)
point3 <- minorep(M=simM3, C=C)

plot(pred ~ Mvec,
     type="n",
     xlim=c(50,75),
     ylim=c(-0.1,1.1),
     ylab="Pr(Minoriy Electoral Success)",
     xlab="M (Raical Margin of Victory)",
     main="",
     mgp=c(2,0.7,0), 
     cex.lab=1.2)
abline(h=c(0, 0.5, 1), lty=2, col="gray60", lwd=1.5)
lines(pred ~ Mvec,
      lwd=4,
      col="maroon",)
text(x=52, y=0.85, labels="C=40", font=2, cex=1.2, col="dimgray")
text(x=56, y=0.75, labels="Pr(minority bloc) = 1", font=2, col="dimgray")
points(x=simM1, y=point1, cex=1.5, col="maroon", pch=16)
points(x=simM2, y=point2, cex=1.5, col="maroon", pch=16)
points(x=simM3, y=point3, cex=1.5, col="maroon", pch=16)
text(x=simM1+7, y=point1+0.02, labels="Pr(crossover) = 0.3", font=2, col="dimgray")
text(x=simM2+6, y=point2, labels="Pr(crossover) = 0.4", font=2, col="dimgray")
text(x=simM3-2, y=point3+0.08, labels="Pr(crossover) = 0.5", font=2, col="dimgray")
title("A", adj=0)


#++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++# 
# (Panel B) STRONG Bloc VOTING
#++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++#
C <- seq(from=1, to=100, by=0.1)
start <- 42
end <- 52

sim1 <- redistrict(coethnic=1, crossover=0)
sim2 <- redistrict(coethnic=1, crossover=0.35)

plot(0, type="n", 
     ylim=c(-0.1,1.1),
     xlim=c(start-5,end+5),
     ylab="Pr(Minority Electoral Success)",
     xlab="C (% of Minority Voters)",
     mgp=c(2,0.7,0),
     cex.lab=1.2)
abline(v=50, lty=2, col="gray60", lwd=1.5)
lines(sim1 ~ C, col="seagreen",lwd=4)
lines(sim2 ~ C, col="maroon",lwd=4)
rect(35, -0.2, start, 1.2, col=alpha("gray80", 0.5), lty=0)
rect(end, -0.2, 65, 1.2, col=alpha("gray80", 0.5), lty=0)
points(x=start, y=sim1[C==start], pch=1, cex=2, col=scales::alpha("seagreen",0.9))
points(x=end, y=sim1[C==end], pch=16, cex=2, col=scales::alpha("seagreen",0.9))
text(x=start, y=sim1[C==start]-0.09, labels=round(sim1[C==start],d=3), col="seagreen")
text(x=end+1, y=sim1[C==end]-0.09, labels=round(sim1[C==end],d=3), col="seagreen")
points(x=start, y=sim2[C==start], pch=1, cex=2, col=scales::alpha("maroon", 0.9))
points(x=end, y=sim2[C==end], pch=16, cex=2, col=scales::alpha("maroon", 0.9))
text(x=start-1, y=sim2[C==start]+0.09, labels=round(sim2[C==start],d=3), col="maroon")
text(x=end, y=sim2[C==end]+0.09, labels=round(sim2[C==end],d=3), col="maroon")
text(x=start+2, y=1.1, labels="Pr(crossover) = 0.35",
     cex=1, col="maroon", font=2)
text(x=start+10, y=-0.09, labels="Pr(crossover) = 0",
     cex=1, col="seagreen", font=2)
arrows(x0=42,x1=45,y0=0,y1=0, length=0.1, lwd=3, col="seagreen")
arrows(x0=44,x1=45,y0=1,y1=1, length=0.1, lwd=3, col="maroon")

 
title("B", adj=0)
 

#++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++# 
# (Panel C) MODERATE Bloc VOTING
#++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++#
 
threshold <- 0.80
C.pr <- 70 # SPECIFIC PERCENTAGE OF CONCERN
sim1 <- redistrict(coethnic=0.9, crossover=0.2)

plot(0, type="n", ylim=c(-0.1,1.1),xlim=c(30,70),
   ylab="Pr(Minority Electoral Success)",xlab="C (% of Minority Voters)",
   mgp=c(2,0.7,0), cex.lab=1.2)
lines(sim1 ~ C, col="maroon",lwd=4)
abline(h=threshold, lty=2, col="dimgray")
text(x=37,y=threshold+0.1, labels="Pre-specified \nThreshold", font=2, col="dimgray")

sw1 = min(C[sim1>=0.8]) # SWEET SPLOT for sim1
rect(sw1, -0.2, C.pr, 1.2, col=alpha("gray80", 0.5), lty=0)
abline(v=C.pr, lwd=2, col="gray60", lty=2)
arrows(x0=sw1, x1=sw1,
       y0=threshold, y1=-0.05, col="maroon", lwd=1, length=0.1, lty=1)
arrows(x0=sw1, x1=sw1,
       y0=0, y1=-0.05, col="maroon", lwd=1, length=0.1)
points(x=sw1, y=threshold, cex=1.5, pch=16, col="maroon")
text(x=sw1+6, y=-0.1, labels=paste0(sw1, "% (Sweet Spot)"), col="maroon", font=2)

text(x=59, y=0.6, labels="Degree of \nPotential Vote Dilution \n(C'-Sweet Spot)",
      cex=1, col="dimgray", font=1)
text(x=64, y=0.2, labels="C'\n(District Plan \nof Interest)",
      cex=1, col="dimgray", font=2)
arrows(x0=65.5, x1=69,
       y0=0.28, y1=0.28, col="dimgray", lwd=1, length=0.1)
title("C", adj=0)


#++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++# 
# (Panel D) Predicting the Number of Minority Officeholders
#++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++#

Mvec <- c(50,40,40,35,70,85)
Cvec <- c(50,40,60,30,50,80)

pred <- minorep(M=Mvec, C=Cvec)
num <- n.minorep(pred)

hist(num,
     ylab="Samples from Monte Carlo Simulations",
     xlab="Predicted Number of Minority Officeholders",
     lwd=1,
     border=F,
     col="gray40",
     main="",
     mgp=c(2,0.7,0),
     cex.lab=1.2)
title("D", adj=0)




dev.off() # Close the Pdf File




#########################################################################################
# END OF THIS R SOURCE FILE
#########################################################################################
